/*
 * Tema 2 ASC
 * 2020 Spring
 */
#include "utils.h"

/*
 * Add your BLAS implementation here
 */

#include <string.h>
#include <cblas.h>

double* my_solver(int N, double *A, double *B) {
	int i, j;
	double *C, *D, *E;

	C = (double *) malloc(N * N * sizeof(double));
	D = (double *) malloc(N * N * sizeof(double));
	E = (double *) malloc(N * N * sizeof(double));

	memcpy(D, B, N * N * sizeof(double));

	cblas_dtrmm(
		CblasRowMajor,	/* 	 CBLAS_LAYOUT 	*/
		CblasRight,		/*    CBLAS_SIDE	*/
		CblasUpper,		/*    CBLAS_UPLO	*/
		CblasTrans,		/*  CBLAS_TRANSPOSE */
		CblasNonUnit,	/*    CBLAS_DIAG	*/
		N,				/* 		  M		    */
		N,				/* 		  N 		*/
		1.0,  			/* 		ALPHA 		*/
		A,				/* 		  A	 		*/
		N,				/* 		 LDA 		*/
		D,				/* 		  B 		*/
		N				/* 		 LDB		*/
	);

	memcpy(E, B, N * N * sizeof(double));

	cblas_dtrmm(
		CblasRowMajor,	/* 	 CBLAS_LAYOUT 	*/
		CblasLeft,		/*    CBLAS_SIDE	*/
		CblasUpper,		/*    CBLAS_UPLO	*/
		CblasNoTrans,	/*  CBLAS_TRANSPOSE */
		CblasNonUnit,	/*    CBLAS_DIAG	*/
		N,				/* 		  M		    */
		N,				/* 		  N 		*/
		1.0,  			/* 		ALPHA 		*/
		A,				/* 		  A	 		*/
		N,				/* 		 LDA 		*/
		E,				/* 		  B 		*/
		N				/* 		 LDB		*/
	);

	cblas_dtrmm(
		CblasRowMajor,	/* 	 CBLAS_LAYOUT 	*/
		CblasLeft,		/*    CBLAS_SIDE	*/
		CblasUpper,		/*    CBLAS_UPLO	*/
		CblasNoTrans,	/*  CBLAS_TRANSPOSE */
		CblasNonUnit,	/*    CBLAS_DIAG	*/
		N,				/* 		  M		    */
		N,				/* 		  N 		*/
		1.0,  			/* 		ALPHA 		*/
		A,				/* 		  A	 		*/
		N,				/* 		 LDA 		*/
		E,				/* 		  B 		*/
		N				/* 		 LDB		*/
	);

	for (i = 0; i < N; i++)
		for (j = 0; j < N; j++)
			C[(i * N) + j] = D[(i * N) + j] + E[(i * N) + j];

	free(D);
	free(E);

	return C;
}
